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The resonance lower-hybrid wave-ion interaction is described by a group (differentiable map) of transfor- 
mations of phase space of the system. All solutions to the map belong to a strange attractor, and chaotic 
motion of the attractor manifests itself in a number of macroscopic effects, such as the energy spectrum 
and particle heating. The applicability of the model to the problem of ion heating by waves at the front of 
collisionless shock as well as ion acceleration by a spectrum of waves is discussed. 
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1. Introduction 

The motion of ions in the lower-hybrid (LH) wave may be the reason for the enhanced particle 
and energy transport, and the interest in ion cyclotron heating is due to the possibility of using 
this plasma heating method in fusion devices (7). The problem is also closely associated with the 
second-order Fermi acceleration of particle at the front of collisionless shocks { 2 - 4 ). There is at 
present some interest in beat-wave acceleration schemes which are able to heat ions directly. So, 
recently Benisti et al. (5) described the nonlinear ion acceleration by a pair of electrostatic waves 
in the LH range. 

First, the problem of ion cyclotron heating has been considered by Karney (6) and discussed 
in (7, Chapter 2). A canonical perturbation theory has been used to examine the transition from 
regular to stochastic motion. Karney showed that the ion gains energy only stochastically, when 
the magnitude of electric field exceeds some threshold value, and using the general notions of 
overlapping resonances and numerical calculations, the threshold value has been estimated. 

The goal of this work is to describe the ion motion in the LH modes under conditions of 
global chaos using an adiabatic approach. We present numerical calculations for particle motion 
in a stochastic regime described in terms of a group of transformations of phase space. The 
dependence of the upper bound of the energy spectrum, heating rate and cross-field diffusion on 
the amplitude of wave field are presented. 

In Section 2, we discuss the equations of motion, and conditions under which chaotic motion 
to occur. In Section 3, we present dynamics of the system as a group of transformations of phase 
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space. In Section 4, we described the strange attractor (SA) of dynamical equations on which the 
motion appeared chaotic. Stochastic diffusion is studied by a Fokker-Plank-Kolmogorov (FPK) 
equation in Section 5. Relying on the results, we consider some possible applications in Section 
6. In Section 7, we give the conclusion of our studies. 


2. Basic equations: the nonlinear resonance wave-particle interaction 

We will discuss dynamics of an ion gyrating in a uniform magnetic field B, and interacting with 
an electrostatic field n(r, 1) of fast magnetohydrodynamic wave at the LH mode, 

© = «LH = (^e^) V2 , (1) 

propagating across the external field. The Hamiltonian of the problem is 

- 0 u(r , t) 

H(r,p,t) = (p-A) 2 + -±-L, (2) 

2m i 

where w; is the ion mass, p the canonical momentum, A the vector potential for the field B and 
f2j(£!i) the ion (electron) gyrofrequency, and we have employed here the system of units in which 
the charge and velocity of light are equal to 1 . 

One chooses a Cartesian system coordinate with the CF-axis directed along 

B, B — (0, 0, B), r — (x, y, z), 

and writes down the vector potential and wave field as 

A = (0, B x , 0), (3) 

u(r, t) — Uq cos( kx — cut). (4) 


Taking into account the axial symmetry of unperturbed system, we introduce the action-angle 
variables (/, 0) carrying through the canonical transformations 

x — r cos 6, y = —mircos sin 6, 

( 21 y\ 2 B (5) 

r — , cu B = £2; = — ■ 


\777CUb / llli 

In the variables, the transformed Hamiltonian (2) is found to be 


H{9, /, t) = H 0 (I) + UoXJnikr) cos(nd - cot), 
Hq — cub/, 


( 6 ) 

(7) 


where Hq is the Hamiltonian of unperturbed system. In deriving Equation (6), we have used the 
Bessel expansion 


cxpiikr cos 6 — cut) = Y, n J n (kr) exp {in6 — cut), n e Z, 


where Z is the set all integers. 
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The equations of motion with the Hamiltonian are 


= U 0 nJ„(- ) sin i/r. 


SH 

0 — — — (9) 

01 

i jf — nO — cot, (10) 

where i/f is the phase of wave, and we take into account only the leading term in Equation (9). 

As is known, the motion of non-autonomous nonlinear dynamics system is mainly determined 
by the behavior of a system near its resonances (7, 8). In first order, the perturbation excites only 
resonance between the frequency u> and the various harmonics &>b, so according to Equations (4), 
(9) and (10) the resonance condition is 

i/f = oo (/) = co — s&>b = 0, s e {s( e Z (1 1) 

co — kv. (12) 

Condition (12) resulting from Equation (11) means that only particles with velocity comparable 
to a phase velocity of wave take part in the resonance wave-particle interaction. 

Again we reveal from Equation (11) that the given system is intrinsically degenerate because 
I is the action variable and Sco(I)/8I = 0 (9). As a consequence, we have to take into account in 
Equation (9) a nonlinear frequency shift acquired by a particle under wave-particle interaction. 
Applying Equation (6) to Equation (9), we have 


«nl(I) = U 0 T, n kJ n (kr ) 


mco^I 


where 7„(-) = 8J n (-)/8kr and relation (5) is used. 

The term takes off the degeneracy, and now the resonance phase space has a typical structure 
of trivial fibering B x S, where B is the base, S the fiber and pair C//, 7 S ) are the coordinates on 
B x S, I e S, \[r (mod) e S, the circle, and each fiber is given by the condition 


s I «b + kU 0 J s (s ) 


~ I i r; i T i I w \ ^ ■ / 

y \TYlG)j$I s J J 

We calculated the distance, 81 = I s+ i — I s , between vth and ( s + l)th resonance fiber using 
Equation (14), which, with s 1, yields 


Uri COnI — mv„ 2 

<*>=■=?, E s = -^— (16) 

Cj S l 

where E s is the particle energy in an exact resonance. Considering Equation (7), the variation of 
I in the vicinity of resonance is given by 

(A/) = sU 0 Js(s) sin i/f, A I = I - I s . (17) 

When one integrates the equation over one period, T = 1 /o> K , we find the width of resonance 
fiber 

^ = s<t>J s (s). (18) 
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It is clear, if the required 


A / SI 
— > — 

Is - I 


( 19 ) 


is valid, that the trivial structure of fibering would be destroyed and the phase space is modified. 
It is obvious that the test is an analog of the Chirikov overlap criterion, which describes the onset 
of stochasticity (10). We will show later that phase space corresponding to chaotic motion is to 
be an SA, which in general is not a topological object. 

At last, substituting Equations (15) and (18) into Equation (19), and using the asymptotic values 
ofBessel functions, J s (s) ~ l\s'^ 3 , J S (-)J^(-) ~ l/s,,y 1 , we obtain the magnitude of the wave 
field 

0 > , ( 20 ) 

s 

at which we should expect the appearance of chaotic phase trajectories. 

Note, because o>nl given by Equation (13) is a decreasing function of /, and expressions (15) 
and (18) have also different dependencies on /, the range of obtainable values of I will be limited. 
This effect will be discussed in detail in Sections 3 and 4. 


3. Dynamics of the system as a recurrent process 


We describe the motion of an ion interacting resonantly with a high-frequency wave. In the 
situation, the characteristic time of wave-particle collision, f c ~ I /w B , is much larger than the 
period of wave oscillations, r w ~ l/&>, 

fc » t w , (21) 

and evolution of the particle phase is a slow process, such that 

\p cot. (22) 


This allows us to treat dynamics of given system in an adiabatic approach, assuming that the 
small parameter s, 


s _ Wb 


(23) 


cu 

serves as a condition for the applicability of the approach. Keeping in mind Equation (7), we carry 
out the following transformations in the RHS of Equation (7). First, we write down 


Uo'EJ n (kr) exp (inO — cot) — UqsJ s (s) sin ip E„^ s exp(im» B f), (24) 

where ip — s6 — cot is the slow variable, whose derivative tends to zero as the system approaches 
an exact resonance. Applying to Equation (24) the Poisson sum formula 

2tt 

£ expO'uftJef) = TT,8(t — nT), T = — , (25) 

coB 

where S(-) in the Dirac delta function, substituting Equations (24) and (25) into Equation (7), we 
arrive at the equation 

E — s 2 ^@E s S(t — nT) sin ip, (26) 

written with the preceding notation (16). 

Additionally, we define E — oo B I, where E is the particle energy, and have used the following 
expression for Bessel asymptotes, J s (s) ~ J r s (s) = J s (s)/s 1 ^. Then one takes into 
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account the particle phase in the exact resonance is equal to 0, and &> NL is given by Equation (13), 
we accomplish just the same procedure with the equation ^ = sO — &>, to obtain 

ir = sco B + Q) (^) s* V8(f - nT) - co. (27) 

Having integrated these equations one by one, we get the closed set of nonlinear difference 
equations 


| I — E n 

fn + 1 = fn 



(mod 2i r) 


(28) 

(29) 


where E n and i j/ n are the variables given at the moment t = nT . 

To study the onset the stochasticity, we linearize Equations (28) and (29) about E„ — E s , 
introducing the new variable £ , 


E n 

l + Sn. 

E s 

Substitution of Equation (30) converts Equations (28) and (29) to the map 


/«\ 2/3 

Cn+I = f n + O Sin l/r„ , 

( 1 \ /CO\ 4 / 3 

y„+i = y n + 7 ( - ®C«+t, (mod 27T) 


4 V©/ 


It should be noted that all maps obtained in this section have the a functional form, 


(30) 


(31) 

(32) 


u„+i — u n + f n+ i = f n + f(n„+i) (mod 2jt), 


and therefore describe dynamics as a recurrent process. Indeed, putting an initial condition, for 
example, uq, at n = 0, at once we get the one-step map, u\ — uq + F(to), i/q = t/t 0 + i jr(u\), 
and so on for all values of n . Solutions of these equations will be studied in the following section. 


4. Strange attractor 


First, we address map (28) and (29). One rewrites these equations in the wave frame. To do this, 
we introduce the new variable n. 


u 


W\V 

m , 

2 E ph 


(33) 


where | V \ is the magnitude of particle velocity. 

In this representation, Equations (28) and (29) go over into the map, g(U, t/r), 


u n+ 1 = u„ + Q sin f n 


(34) 



0K-ir 1/2 


1 

\jf n+ i = \lr n + - 


2/3 


(mod 2: r ) , 


(35) 
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where the control parameter Q is given by 


Q = 





It is obvious that (t/r, u) are the variables on a 2D smooth manifold M, and 


( 36 ) 


g:M 


M 


is the family of maps for all values of u,n e Z. 
By virtue of the following properties of g. 


g n+1 = g n g l ,g n = (g 1 r, 


(37) 


where g 1 is the first-step map, g is the group of transformation of the phase space M and, 
consequently, the pair (M, g) is a dynamical system given by initial conditions. 

To prove that the system demonstrates a chaotic motion, we need to show that all solutions 
irrespective of initial conditions belong to an SA. Then the eigenvalues 7, | and 7,2 of the Jacobian 
matrix, 7: 


9 (^h«+i , u n +i ) 
9('h/ ! , w n ) 


(38) 


are to be found. Denote by 

det 7 = 7.2, tr 7 = A.i -(- 7.2 


(39) 


the determinant and trace of this matrix. 

Applying Equations (34) and (35) to Equation (38) yields 


det 7=1. 


(40) 


Therefore, this g is the measure-preserving map, and the pair (i//, u) is the canonical pair of 
variables. 

It is known (9) that the requirement | tr 7 1 — 1 > 2 defines topological equivalence of hyperbolic 
sets, and the relation 

|tr 7| — 1 = 2 (41) 

is the condition of topological modification of phase space. Thus, 


|tr 7 1 = 2 + 




this condition allows us to calculate the upper bound of {«} 


sup{;r} = \ub\ = 




(42) 


(43) 


The formula predicts the Q 4 ^ - dependence of a j, on the control parameter Q. 

We have numerically integrated Equations (34) and (35) for several different values of Q. The 
initial conditions were chosen in a random fashion and corresponded to the region of small values 
of (i/r, u). Figure 1 shows some of our results computed for one trajectory in the i jr — u phase 
space. The figure shows that the boundary of {u} is well approximated by formula (43). 
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Figure 1. The phase space of the map. One single trajectory of length 10 4 for different values of Q. (a) Q = 16 and 
(b) Q = 32. 


The existence of upper bound implies that the given manifold is a compact space, on which we 
can determine the eigenvalues 7i and 7, 2 , the Kolmogorov entropy, K, being a kneading invariant, 


K = lnA.i, 


and the fractal dimension of the M, 


d{=l- 


ln 7.1 
In 7.2 ' 


(44) 


(45) 


By calculating, we found K — ln[(3 + V5)/2]; therefore, the rate of a loss of information 
is positive. Then we compute the fractal measure df = 1 , implying that the phase curve whose 
topological measure is equal to one evenly fills all obtainable phase space. By virtue of these 
invariants, the pair (M, g) is an SA tightly embedded in the phase space at n — »■ 00 . 

We have found in Section 3 that transition of a given system to chaotic motion can be described 
by the set of difference equations ((31) and (32)). 

It is clear, that by rescaling f (l/4)(&>/a>„) 4 / 3 <t>§, these equations are reduced to the map: 


1 / a> \ 2 

£«+ 1 = Hn + T — <f> sm'I', 

4 \(o B J 

T'n+i = + ^ n+ \ (mod oj). 


(46) 

(47) 


having the form of a standard map (7, 8). We will again use topological methods to study the 
behavior of the map. First, we compute the trace of Jacobian of the map. 


tri = 



When this expression satisfies Equation (41), we find 


4>c 


Uc _ 

T ph (O 


(48) 


( 49 ) 


the minimal value of the wave field required for the onset of stochastic motion. 
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Now substituting Equation (49) into Equation (43), we write down the simple expression for 
the upper boundary of spectrum, 


\Ub\ — 


E 

E p h 



(50) 


Let us consider the dynamics of ions in an MH wave at the LH frequency, &>lh- In this 
case, we must set everywhere in our formulas u>b — £2;, &> = &>lh = where £2; £2 e is 

the gyrofrequency of an ion (electron). 

Finally, Equation (49) takes the form 


U c = 2 



E ph — 


111 i V 


2 

ph 


2 ’ 


(51) 


and as it appears from Equations (50) and (51), the feasible extent of particle heating is 



(52) 


5. Heating of ions by the LH wave 


Introduce the distribution function f{E,t) on the stochastic set, the evolution of which obeys the 
FPK equation. 


Bf(E,t) = 1 3_ v df(E,t) 
dt 2 8E E 8E 

Here 'D v is the coefficient of diffusion in energy space given by 


(53) 


2? £ = 


(E n+ 1 - E„) 2 
T 


(54) 


where (•} is the operator of phase average, E n+ \ — E n is given by the original map (34) and (35) 
and T = 2it / £2; is the ion gyroperiod. 

Making use of the map in Equation (54), we get 


2 ? £ = 2 




(55) 


The function f(E, t) is positive on the SA, and beyond the SA, both the function and its derivative 
are equal to zero. Therefore, the function, or rather the probability density, obeys the norm 



J)dE = 1. 


Under the conditions, the steady-state solution to Equation (53) is 


f(E) = 


1 

Eb E ph 


and the characteristic time for establishing the distribution is found to be 





T. 


(56) 


(57) 


(58) 


We have employed in Equation (58) the relations (50) and (55). 
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The / ( E) as expected does describe a uniform distribution on the S A and determines the energy 
spectrum of a particle. One points out the feature of the distribution near the transition to chaos. 
If <t> -> <t> c , then Eb E p h and f(E) -> S(E b — E p \ ,), S is the Dirac delta function, reflecting 
the bifurcation transition to chaos. Like this, <J> C is the bifurcation parameter. 

The distribution permits calculating the means 


r E » 1 

(E) — I f(E)EdE = , 

J Eph 2 (E b + E ph y 

(E 2 ) = X -{E 2 + E b E ph +E 2 ph , ), 
and the relative level of fluctuations. 


- ( E > 2 ) 

fW) 


(59) 

(60) 


( 61 ) 


for $ <t> c . 

The distribution f(E) as well as all means through the very high level of fluctuations are strong, 
stable due to the global stability of the SA. 

The most interesting feature of the system is the dependence of t& on the magnitude of wave 
field. As it follows from Equation (58), oc (<t>(<J> c )) 2,/3 , while this dependence is typically a 
decreasing function of O. This feature is conditioned by the properties of map itself having an 
explicit dependence on the wave field in the i//- -equation. This property is irreducible due to the 
degeneracy of original linear problem. 

Let us discuss how the system approaches the SA. Introducing the second moment, {E 2 ) — 
f E 2 f(E, t)dE, we integrate Equation (53) to find with the help of Equations (55) and (56), 


d(£ 2 } 

~^T 


= V E , 


(62) 


that is, the particle heating is realized by a Brownian process, Ea^/1, and the heating rate is 
approximately given by 


d£ 

df 


= T~ l 



1/3 

<b(<J> c ) 2 £ ph 



(63) 


Note that the diffusion in energy space is accompanied by the space diffusion across an external 
magnetic field. Indeed, gyroradius r and the particle energy E are directly connected by the relation 

r 2 = 2 E/m£l 2 . (64) 


Utilizing Equation (63) in Equation (64), we have the cross-field coefficient 


d(r 2 ) 

ck 



. m e 

Z>o= — 

ntj 


1/3 


rl ph 
mCliir 


6. Application 

First we address the problem of particle acceleration by waves at the front of perpendicular 
collisionless shock. The shock front is sufficiently steep, and the electrostatic potential is localized 
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at the front, the width of which is typically of the order of an inertial length. In shock surfing, 
the particle accelerates along the shock front under the action of the convective electric field of 
the plasma flow to v ~ Vf in a time r ~ where Vf is the front speed. Then an ion with 

v ~ Vf ~ to/k can be trapped by wave field at the bottom of the potential. Remaining close to the 
potential bottom, the ion is accelerated by the electric field of wave, E w , as ever the condition 

E w > v p B (65) 


is fulfilled (8). 

Again we have identified above the test for electric field of MH waves, which is necessary for 
the acceleration to occur. Thus, using E m = k<f>, Equation (49) results in ( E w ) c = v v B, which 
is similar to Equation (65). Furthermore, we have established that in wave fields satisfying the 
condition, ion motion becomes chaotic and leads to stochastic heating of particles. In accordance 
with Equations (52) and (58), particles are accelerated very rapidly as particles are heated to /7 h 
in about 100 gyroperiods. 

Note, particle simulations (77) indicate that electron acceleration by the upper-hybrid wave can 
occur in this regime. 

At last, it should be noted that the effects related to a finite age and size of the shock and the 
particle loss because scattering can modify the energy spectrum, while the heating rate as a rule 
remains the same (72). 

Ion heating by two electrostatic waves in the LH range was proposed by Benisti et al. (5). The 
authors demonstrated that ions with an arbitrary low initial velocity can be accelerated through a 
nonlinear interaction with a pair of waves that obey beating criterions 


(o 2 -(o i = (o B . 


(66) 


Because of the lack of a threshold for the initial ion velocity, this acceleration scheme could be 
promising to many applications, such as plasma heating in fusion devices and spacecraft plasma 
propulsion. 

Afterwards, a numerical exploration of this mechanism revealed that the test (66) is necessary 
but not sufficient, and a second-order perturbation analysis was carried out to define the domains 
of allowed and forbidden accelerations ( 13). They led also to a conclusion that, despite from 
restrictions, an ion with arbitrary low velocity may benefit from this mechanism. 

Following the authors, we introduce these ideas by considering their model by the method 
mentioned above. 

First we determine the wave field of two waves as 


U ( x , t) — 7/o(cos(A'iv — w\t) + cos (k 2 x — a> 2 1)) 


2 

C 0\ + (02 


— Uq cos(kx — cot ) cos (kx — vt ) , (67) 


k — 


k\ + k 2 


v — 


CO 2 — CO l 


K 2 - Kl 


( 68 ) 


According to (9), we set k — 0, which corresponds to the required u g /t> p )$> 1, where v g is the 
group velocity. Then we introduce the action-angle variables to write equation for 7, 


7 = Uq cos vt nJ„{-) sin i fr, \[r = n6 — cot , 


(69) 


and an appropriate equation for the variable 6 . 
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Now, in view of Equations (67) and (69), resonance conditions satisfy the equations 

• CO l + CO 2 

co — sO = 0, co = kv, co = = &>lh- 


(70) 


As before, we assume that the conditions for adiabatic approach i// < < i jrw, I < < col are met, 
and write down Equation (69) at first as 


• x — ^ 27T 

I = sUq cos vtJ s (s)T sin xfr > S(t — nT ), T — — . 

CO B 


(71) 


then integrating it over one time period reduces it to the difference equation, 

/„+ 1 = /„ + sUoJ s (s)T cos vT sin x[r n . (72) 


Unlike Equation (28), there occurs the term cos vT due to beat oscillations, the mean value of 
which is equal to zero in the space of parameters. 

At last, under the conditions of parametric resonance (9), vT — tc and the condition 


a >2 — co\ = cob 


(73) 


is valid. Both Equations (70) and (73) serve as the sufficient and necessary condition for efficient 
interaction, and Equation (72) and an appropriate equation for the particle phase reduce to the set, 
which is tantamount to (28) and (29) entirely. 

Proceed to the contrary case, when v g /v v << 1. For the case, we have in Equation (69) the 
term cosU'i’„'F) instead of cos v T . Taking the dispersion relation 


co = «lh- 


ko= — , 
c 


k 2 
K o 


2 k 3 


Va = 


k » ko, 




2 k 3 


(74) 

(75) 


for a given branch of waves (7), whence, it follows M 2 — M\ — v g (k 2 — k\ ) for m >> &>2 — ft>i. 
Substituting it in the beat oscillation term, cos(v g T(k 2 — Aq)/2), we again arrive at the condition 
for parametric resonance coinciding with Equation (73). So we found that the acceleration scheme 
is identical to that given above, and therefore, it is realized, if and only if condition (49) holds. 
Thus, it has no advantage over the traditional one. 


7. Conclusion 

Using the adiabatic approach, the equations of ion motion in electrostatic field of an LH wave 
reduce to a map. It is shown that all solutions of the map belong to an SA, and the means 
(observables) on the attractor are stable, irrelevant of initial conditions. The solutions have revealed 
that the condition for the onset of global stochasticity is E c = v v B, in close agreement with the 
numerical value. The estimate is larger about (mi/m c ) 1 ^ 6 times than that obtained by Karney (6). 
The upper bound of the SA can be put in correspondence with the upper boundary of energy 
spectrum whose value depends only on the amplitude of wave. Chaotic motion gives rise to 
diffusion in energy and leads to establishing the energy spectrum, and the timescales of the 
process is of the order of tens of gyroperiods. The results have been applied to a number of 
problems, including Fermi acceleration by waves at the front of shocks, and ion cyclotron heating 
in the beat-wave accelerator. 
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